! This module includes a "hill climbing" search routine to find peak of (conditional) value functions
! and a "search" routine that locates the relevant (asset) node in a grid
! and a cubic interpolation function (takes 4 points)

module searchnew

use params
use states 

implicit none

contains

!*******************************************************
  real(prec) function golden_cons(ax,bx,cx,f,xmin,kappaamt,leisure,famcare,insured,cp,ck,vp,vk,cash,nodek,medicaid)
  
  real(prec) f
  real(prec) xmin
  real(prec) kappaamt
  real(prec) leisure
  real(prec) cash

  integer famcare
  integer insured
  integer medicaid

  integer nodek

  integer count

  real(prec) ax,bx,cx
  real(prec) R,C, tol

  real(prec) cp,ck,vp,vk

  EXTERNAL f

  PARAMETER (R=0.61803399d0, C=1.d0-R, tol=1.d-3)

  real(prec) f1,f2,x0,x1,x2,x3

  x0=ax
  x3=cx

  count=0

  if(abs(cx-bx).gt.abs(bx-ax))then
    x1=bx
    x2=bx+C*(cx-bx)
  else
    x2=bx
    x1=bx-C*(bx-ax)
  endif

  f1=f(x1,kappaamt,leisure,famcare,insured,cp,ck,vp,vk,cash,nodek,medicaid)
  f2=f(x2,kappaamt,leisure,famcare,insured,cp,ck,vp,vk,cash,nodek,medicaid)
  count=count+1

1     if(abs(x3-x0).gt.tol*(abs(x1)+abs(x2)))then
        if(f2.lt.f1)then
          x0=x1
          x1=x2
          x2=R*x1+C*x3

          f1=f2
          f2=f(x2,kappaamt,leisure,famcare,insured,cp,ck,vp,vk,cash,nodek,medicaid)
        else
          x3=x2
          x2=x1
          x1=R*x2+C*x0
          f2=f1
          f1=f(x1,kappaamt,leisure,famcare,insured,cp,ck,vp,vk,cash,nodek,medicaid)
        endif

        count=count+1

        goto 1

      endif

  if(f1.lt.f2)then
    golden_cons=f1
    xmin=x1
  else
    golden_cons=f2
    xmin=x2
  endif

  return

  end function golden_cons
! *******************************************************

!*******************************************************
  real(prec) function golden_cons_deadparent(ax,bx,cx,f,xmin,leisure,cash,nodek)
  
  real(prec) f
  real(prec) xmin
  real(prec) leisure
  real(prec) cash

  integer nodek

  integer count

  real(prec) ax,bx,cx
  real(prec) R,C, tol

  real(prec) cp,ck,vp,vk

  EXTERNAL f

  PARAMETER (R=0.61803399d0, C=1.d0-R, tol=1.d-3)

  real(prec) f1,f2,x0,x1,x2,x3

  x0=ax
  x3=cx

  count=0

  if(abs(cx-bx).gt.abs(bx-ax))then
    x1=bx
    x2=bx+C*(cx-bx)
  else
    x2=bx
    x1=bx-C*(bx-ax)
  endif

  f1=f(x1,leisure,cash,nodek)
  f2=f(x2,leisure,cash,nodek)
  count=count+1

1     if(abs(x3-x0).gt.tol*(abs(x1)+abs(x2)))then
        if(f2.lt.f1)then
          x0=x1
          x1=x2
          x2=R*x1+C*x3

          f1=f2
          f2=f(x2,leisure,cash,nodek)
        else
          x3=x2
          x2=x1
          x1=R*x2+C*x0
          f2=f1
          f1=f(x1,leisure,cash,nodek)
        endif

        count=count+1

        goto 1

      endif

  if(f1.lt.f2)then
    golden_cons_deadparent=f1
    xmin=x1
  else
    golden_cons_deadparent=f2
    xmin=x2
  endif

  return

  end function golden_cons_deadparent
! *******************************************************


! *******************************************************
subroutine find_ass(pertype,fass,fage,fmina, wgt1,wgt2)

    integer pertype, fage, fmina, jl, jm, ju , stat
    real(prec) fass, a1,a2, wgt1,wgt2
    real(prec) assmat(na,nt)

    if (pertype.eq.1) then
      assmat=assetp
    else
      assmat=assetk
    endif

    jl=0
    ju=na+1
812  if (ju-jl.gt.1) then
      jm = (ju+jl)/2
      if (fass.gt.assmat(jm,fage) ) then
          jl = jm
      else
          ju = jm
      endif
      goto 812
    endif
    fmina = jl

    !find how far the asset amount lies between two nodes
    if (fmina.eq.na) then         !fix boundary issues
      fmina = na-1
    elseif (fmina.eq.0) then
      fmina = 1
    endif
    a1 = assmat(fmina,fage)      !get asset value at lower node
    a2 = assmat(fmina+1,fage)    !get asset value at higher node

    wgt1=(fass-a2)/(a1-a2)        !get weight on lower node
    wgt2=(fass-a1)/(a2-a1)        !get weight on higher node

end subroutine find_ass
! *******************************************************




! *******************************************************
subroutine find_valpwgt(val,ins,nodek)

  integer m
  real(prec) val
  integer ins
  integer nodek

  do m=mnode,nm

    !if get all the way through to nm, cant do any better so just assign all weight to parent (m=nm)
    if (valuep_married(age,ins,bignode,m,apnode,nodek)>=val .or. m==nm) then

      mu_padj(age,ins,bignode,apnode,nodek)=paretoweight(m)
      if (valuep_married(age,ins,bignode,m,apnode,nodek)==valuepNT0_married(age,ins,bignode,m,apnode,nodek)) then
        kappa_padj(age,ins,bignode,apnode,nodek)=kappaNT0_married(age,ins,bignode,m,apnode,nodek)
        valuep_padj(age,ins,bignode,apnode,nodek)=valuepNT0_married(age,ins,bignode,m,apnode,nodek)
        valuek_padj(age,ins,bignode,apnode,nodek)=valuekNT0_married(age,ins,bignode,m,apnode,nodek)
        conexpp_padj(age,ins,bignode,apnode,nodek)=conexppNT0_married(age,ins,bignode,m,apnode,nodek)
        conexpk_padj(age,ins,bignode,apnode,nodek)=conexpkNT0_married(age,ins,bignode,m,apnode,nodek)
        LFP_padj(age,ins,bignode,apnode,nodek)=1
        famcare_padj(age,ins,bignode,apnode,nodek)=0        
      elseif (valuep_married(age,ins,bignode,m,apnode,nodek)==valuepNT100_married(age,ins,bignode,m,apnode,nodek)) then
        kappa_padj(age,ins,bignode,apnode,nodek)=kappaNT100_married(age,ins,bignode,m,apnode,nodek)
        valuep_padj(age,ins,bignode,apnode,nodek)=valuepNT100_married(age,ins,bignode,m,apnode,nodek)
        valuek_padj(age,ins,bignode,apnode,nodek)=valuekNT100_married(age,ins,bignode,m,apnode,nodek)
        conexpp_padj(age,ins,bignode,apnode,nodek)=conexppNT100_married(age,ins,bignode,m,apnode,nodek)
        conexpk_padj(age,ins,bignode,apnode,nodek)=conexpkNT100_married(age,ins,bignode,m,apnode,nodek)
        LFP_padj(age,ins,bignode,apnode,nodek)=1
        famcare_padj(age,ins,bignode,apnode,nodek)=1
      elseif (valuep_married(age,ins,bignode,m,apnode,nodek)==valuepPT0_married(age,ins,bignode,m,apnode,nodek)) then
        kappa_padj(age,ins,bignode,apnode,nodek)=kappaPT0_married(age,ins,bignode,m,apnode,nodek)
        valuep_padj(age,ins,bignode,apnode,nodek)=valuepPT0_married(age,ins,bignode,m,apnode,nodek)
        valuek_padj(age,ins,bignode,apnode,nodek)=valuekPT0_married(age,ins,bignode,m,apnode,nodek)
        conexpp_padj(age,ins,bignode,apnode,nodek)=conexppPT0_married(age,ins,bignode,m,apnode,nodek)
        conexpk_padj(age,ins,bignode,apnode,nodek)=conexpkPT0_married(age,ins,bignode,m,apnode,nodek)
        LFP_padj(age,ins,bignode,apnode,nodek)=2
        famcare_padj(age,ins,bignode,apnode,nodek)=0        
      elseif (valuep_married(age,ins,bignode,m,apnode,nodek)==valuepPT100_married(age,ins,bignode,m,apnode,nodek)) then
        kappa_padj(age,ins,bignode,apnode,nodek)=kappaPT100_married(age,ins,bignode,m,apnode,nodek)
        valuep_padj(age,ins,bignode,apnode,nodek)=valuepPT100_married(age,ins,bignode,m,apnode,nodek)
        valuek_padj(age,ins,bignode,apnode,nodek)=valuekPT100_married(age,ins,bignode,m,apnode,nodek)
        conexpp_padj(age,ins,bignode,apnode,nodek)=conexppPT100_married(age,ins,bignode,m,apnode,nodek)
        conexpk_padj(age,ins,bignode,apnode,nodek)=conexpkPT100_married(age,ins,bignode,m,apnode,nodek)
        LFP_padj(age,ins,bignode,apnode,nodek)=2
        famcare_padj(age,ins,bignode,apnode,nodek)=1
      elseif (valuep_married(age,ins,bignode,m,apnode,nodek)==valuepFT0_married(age,ins,bignode,m,apnode,nodek)) then
        kappa_padj(age,ins,bignode,apnode,nodek)=kappaFT0_married(age,ins,bignode,m,apnode,nodek)
        valuep_padj(age,ins,bignode,apnode,nodek)=valuepFT0_married(age,ins,bignode,m,apnode,nodek)
        valuek_padj(age,ins,bignode,apnode,nodek)=valuekFT0_married(age,ins,bignode,m,apnode,nodek)
        conexpp_padj(age,ins,bignode,apnode,nodek)=conexppFT0_married(age,ins,bignode,m,apnode,nodek)
        conexpk_padj(age,ins,bignode,apnode,nodek)=conexpkFT0_married(age,ins,bignode,m,apnode,nodek)
        LFP_padj(age,ins,bignode,apnode,nodek)=3
        famcare_padj(age,ins,bignode,apnode,nodek)=0        
      elseif (valuep_married(age,ins,bignode,m,apnode,nodek)==valuepFT100_married(age,ins,bignode,m,apnode,nodek)) then
        kappa_padj(age,ins,bignode,apnode,nodek)=kappaFT100_married(age,ins,bignode,m,apnode,nodek)
        valuep_padj(age,ins,bignode,apnode,nodek)=valuepFT100_married(age,ins,bignode,m,apnode,nodek)
        valuek_padj(age,ins,bignode,apnode,nodek)=valuekFT100_married(age,ins,bignode,m,apnode,nodek)
        conexpp_padj(age,ins,bignode,apnode,nodek)=conexppFT100_married(age,ins,bignode,m,apnode,nodek)
        conexpk_padj(age,ins,bignode,apnode,nodek)=conexpkFT100_married(age,ins,bignode,m,apnode,nodek)
        LFP_padj(age,ins,bignode,apnode,nodek)=3
        famcare_padj(age,ins,bignode,apnode,nodek)=1
      endif

      goto 1234
    endif

  enddo

1234 continue
end subroutine find_valpwgt
! *******************************************************


! *******************************************************
! *******************************************************
subroutine find_valkwgt(val,ins,nodek)

  integer m
  real(prec) val
  integer ins
  integer nodek

  do m=mnode,1,-1

    !if get all the way through to 1, cant do any better so just assign all weight to kid (m=1)
    if (valuek_married(age,ins,bignode,m,apnode,nodek)>=val .or. m==1) then

      mu_kadj(age,ins,bignode,apnode,nodek)=paretoweight(m)
      if (valuek_married(age,ins,bignode,m,apnode,nodek)==valuekNT0_married(age,ins,bignode,m,apnode,nodek)) then
        kappa_kadj(age,ins,bignode,apnode,nodek)=kappaNT0_married(age,ins,bignode,m,apnode,nodek)
        valuep_kadj(age,ins,bignode,apnode,nodek)=valuepNT0_married(age,ins,bignode,m,apnode,nodek)
        valuek_kadj(age,ins,bignode,apnode,nodek)=valuekNT0_married(age,ins,bignode,m,apnode,nodek)
        conexpp_kadj(age,ins,bignode,apnode,nodek)=conexppNT0_married(age,ins,bignode,m,apnode,nodek)
        conexpk_kadj(age,ins,bignode,apnode,nodek)=conexpkNT0_married(age,ins,bignode,m,apnode,nodek)
        LFP_kadj(age,ins,bignode,apnode,nodek)=1
        famcare_kadj(age,ins,bignode,apnode,nodek)=0        
      elseif (valuek_married(age,ins,bignode,m,apnode,nodek)==valuekNT100_married(age,ins,bignode,m,apnode,nodek)) then
        kappa_kadj(age,ins,bignode,apnode,nodek)=kappaNT100_married(age,ins,bignode,m,apnode,nodek)
        valuep_kadj(age,ins,bignode,apnode,nodek)=valuepNT100_married(age,ins,bignode,m,apnode,nodek)
        valuek_kadj(age,ins,bignode,apnode,nodek)=valuekNT100_married(age,ins,bignode,m,apnode,nodek)
        conexpp_kadj(age,ins,bignode,apnode,nodek)=conexppNT100_married(age,ins,bignode,m,apnode,nodek)
        conexpk_kadj(age,ins,bignode,apnode,nodek)=conexpkNT100_married(age,ins,bignode,m,apnode,nodek)
        LFP_kadj(age,ins,bignode,apnode,nodek)=1
        famcare_kadj(age,ins,bignode,apnode,nodek)=1
      elseif (valuek_married(age,ins,bignode,m,apnode,nodek)==valuekPT0_married(age,ins,bignode,m,apnode,nodek)) then
        kappa_kadj(age,ins,bignode,apnode,nodek)=kappaPT0_married(age,ins,bignode,m,apnode,nodek)
        valuep_kadj(age,ins,bignode,apnode,nodek)=valuepPT0_married(age,ins,bignode,m,apnode,nodek)
        valuek_kadj(age,ins,bignode,apnode,nodek)=valuekPT0_married(age,ins,bignode,m,apnode,nodek)
        conexpp_kadj(age,ins,bignode,apnode,nodek)=conexppPT0_married(age,ins,bignode,m,apnode,nodek)
        conexpk_kadj(age,ins,bignode,apnode,nodek)=conexpkPT0_married(age,ins,bignode,m,apnode,nodek)
        LFP_kadj(age,ins,bignode,apnode,nodek)=2
        famcare_kadj(age,ins,bignode,apnode,nodek)=0        
      elseif (valuek_married(age,ins,bignode,m,apnode,nodek)==valuekPT100_married(age,ins,bignode,m,apnode,nodek)) then
        kappa_kadj(age,ins,bignode,apnode,nodek)=kappaPT100_married(age,ins,bignode,m,apnode,nodek)
        valuep_kadj(age,ins,bignode,apnode,nodek)=valuepPT100_married(age,ins,bignode,m,apnode,nodek)
        valuek_kadj(age,ins,bignode,apnode,nodek)=valuekPT100_married(age,ins,bignode,m,apnode,nodek)
        conexpp_kadj(age,ins,bignode,apnode,nodek)=conexppPT100_married(age,ins,bignode,m,apnode,nodek)
        conexpk_kadj(age,ins,bignode,apnode,nodek)=conexpkPT100_married(age,ins,bignode,m,apnode,nodek)
        LFP_kadj(age,ins,bignode,apnode,nodek)=2
        famcare_kadj(age,ins,bignode,apnode,nodek)=1
      elseif (valuek_married(age,ins,bignode,m,apnode,nodek)==valuekFT0_married(age,ins,bignode,m,apnode,nodek)) then
        kappa_kadj(age,ins,bignode,apnode,nodek)=kappaFT0_married(age,ins,bignode,m,apnode,nodek)
        valuep_kadj(age,ins,bignode,apnode,nodek)=valuepFT0_married(age,ins,bignode,m,apnode,nodek)
        valuek_kadj(age,ins,bignode,apnode,nodek)=valuekFT0_married(age,ins,bignode,m,apnode,nodek)
        conexpp_kadj(age,ins,bignode,apnode,nodek)=conexppFT0_married(age,ins,bignode,m,apnode,nodek)
        conexpk_kadj(age,ins,bignode,apnode,nodek)=conexpkFT0_married(age,ins,bignode,m,apnode,nodek)
        LFP_kadj(age,ins,bignode,apnode,nodek)=3
        famcare_kadj(age,ins,bignode,apnode,nodek)=0        
      elseif (valuek_married(age,ins,bignode,m,apnode,nodek)==valuekFT100_married(age,ins,bignode,m,apnode,nodek)) then
        kappa_kadj(age,ins,bignode,apnode,nodek)=kappaFT100_married(age,ins,bignode,m,apnode,nodek)
        valuep_kadj(age,ins,bignode,apnode,nodek)=valuepFT100_married(age,ins,bignode,m,apnode,nodek)
        valuek_kadj(age,ins,bignode,apnode,nodek)=valuekFT100_married(age,ins,bignode,m,apnode,nodek)
        conexpp_kadj(age,ins,bignode,apnode,nodek)=conexppFT100_married(age,ins,bignode,m,apnode,nodek)
        conexpk_kadj(age,ins,bignode,apnode,nodek)=conexpkFT100_married(age,ins,bignode,m,apnode,nodek)
        LFP_kadj(age,ins,bignode,apnode,nodek)=3
        famcare_kadj(age,ins,bignode,apnode,nodek)=1
      endif

      goto 1235
    endif

  enddo

1235 continue

end subroutine find_valkwgt
! *******************************************************

! *******************************************************
subroutine find_pwgt(fwgt,fminp, wgt1,wgt2)

    integer fminp
    integer jl, jm, ju
    real(prec) fwgt
    real(prec) p1, p2
    real(prec) wgt1, wgt2

          jl=0
          ju=nm+1
812       if (ju-jl.gt.1) then
            jm = (ju+jl)/2
            if (fwgt.gt.paretoweight(jm)) then
                jl = jm
            else
                ju = jm
            endif
            goto 812
          endif
          fminp = jl

          !find how far the pwgt lies between two nodes
          if (fminp.eq.nm) then         !fix boundary issues
            fminp = nm-1
          elseif (fminp.eq.0) then
            fminp = 1
          endif
          p1 = paretoweight(fminp)    !get pareto weight value at lower node
          p2 = paretoweight(fminp+1)  !get pareto weight value at higher node

          wgt1=(fwgt-p2)/(p1-p2)        !get weight on lower node
          wgt2=(fwgt-p1)/(p2-p1)        !get weight on higher node

end subroutine find_pwgt
! *******************************************************

! *******************************************************
real (prec) function interp_deadparent(iage,inodeto,imina,x,func)

    integer iage
    integer imina,node1,node2,inodeto
    real (prec) x,y1,y2,iwgt1,iwgt2
    real (prec) func(nt,nonodesk)

    !LINEAR INTERPOLATION
    node1=imina
    node2=imina+1

    y1 = func(iage,inodeto+node1)
    y2 = func(iage,inodeto+node2)
    iwgt1 = ((x-node2))/((node1-node2))
    iwgt2 = ((x-node1))/((node2-node1))
    interp_deadparent=iwgt1*y1+iwgt2*y2


end function interp_deadparent
! *******************************************************

! *******************************************************
real (prec) function interp_val(minnode,midnode,func,ins,nodek)

    integer minnode
    integer ins,ins_lb,ins_ub
    real(prec) midnode
    real(prec) func(nt,ni,nbignodes,nm,na,nonodesk)

    integer node1,node2
    real(prec) y1,y2
    real(prec) iwgt1,iwgt2

    integer nodek

    !LINEAR INTERPOLATION
    node1=minnode
    node2=minnode+1

    if (age .eq. 1) then
      ins_lb=insurance_married(bignode,node1,apnode,nodek)
      ins_ub=insurance_married(bignode,node2,apnode,nodek)
      y1 = func(age,ins_lb,bignode,node1,apnode,nodek)
      y2 = func(age,ins_ub,bignode,node2,apnode,nodek)
      iwgt1 = ((midnode-node2))/((node1-node2))
      iwgt2 = ((midnode-node1))/((node2-node1))
      interp_val=iwgt1*y1+iwgt2*y2
    else
      y1 = func(age,ins,bignode,node1,apnode,nodek)
      y2 = func(age,ins,bignode,node2,apnode,nodek)
      iwgt1 = ((midnode-node2))/((node1-node2))
      iwgt2 = ((midnode-node1))/((node2-node1))
      interp_val=iwgt1*y1+iwgt2*y2
    endif

end function interp_val
! *******************************************************

! *******************************************************
real(prec) function bilinear(x, p1, p2, values) result(r)
  ! Returns bilinear interpolation of a point 'x' using uniform data 'values'
  ! specified in a box determined by 'p1' and 'p2' points.
  real(prec), intent(in) :: x(:) ! The 2D coordinates of a point to interpolate
  real(prec), intent(in) :: p1(:) ! The lower left front corner is p1(2)
  real(prec), intent(in) :: p2(:) ! The upper right back corner is p2(2)
  real(prec), intent(in) :: values(:, :) ! Values on a uniform 2D grid
  real(prec) :: x0(2)
  integer :: ij(2), nelem(2)
  ! Number of elements in each direction:
  nelem = shape(values) - 1
  ! Transform to box [0,1] x [0,1] * nelem, that is, make the
  ! element length exactly 1 in each direction
  x0 = nelem * (x - p1) / (p2 - p1)
  ij = int(x0)+1 ! indices of the nearest vertex
  where (ij > nelem) ij = nelem
  where (ij < 1) ij = 1
  r = bilinear_unitbox(x0-ij+1, values(ij(1):ij(1)+1, ij(2):ij(2)+1))
end function

real(prec) function bilinear_unitbox(xy, V) result(Vxy)
  ! Returns bilinear interpolation of a point (x, y) in a unit box
  real(prec), intent(in) :: xy(2)  ! The point in a unit box [0,1] x [0,1]
  ! V(i, j) ... value at vertex i, j of the box, where i,j = 1,2:
  real(prec), intent(in) :: V(:, :)
  real(prec) :: x, y
  x = xy(1)
  y = xy(2)
  Vxy = &
      V(1, 1) * (1-x)*(1-y) + &
      V(2, 1) *   x  *(1-y) + &
      V(1, 2) * (1-x)*  y   + &
      V(2, 2) *   x  *  y
end function
! *******************************************************

! *******************************************************
real(prec) function trilinear(x, p1, p2, values) result(r)

    ! Returns trilinear interpolation of a point 'x' using uniform data 'values'
    ! specified in a box determined by 'p1' and 'p2' points.
    real(prec), intent(in) :: x(:) ! The 3D coordinates of a point to interpolate
    real(prec), intent(in) :: p1(:) ! The lower left front corner is p1(3)
    real(prec), intent(in) :: p2(:) ! The upper right back corner is p2(3)
    real(prec), intent(in) :: values(:, :, :) ! Values on a uniform 3D grid
    real(prec) :: x0(3)
    integer :: ijk(3), nelem(3)
    ! Number of elements in each direction:
    nelem = shape(values) - 1
!     print*,'nelem',nelem
    ! Transform to box [0,1] x [0,1] x [0,1] * nelem, that is, make the
    ! element length exactly 1 in each direction
    x0 = nelem * (x - p1) / (p2 - p1)
    ijk = int(x0)+1 ! indices of the nearest vertex
    where (ijk > nelem) ijk = nelem
    where (ijk < 1) ijk = 1
    r = trilinear_unitbox(x0-ijk+1, values(ijk(1):ijk(1)+1, ijk(2):ijk(2)+1, &
            ijk(3):ijk(3)+1))
end function

real(prec) function trilinear_unitbox(xyz, V) result(Vxyz)
    ! Returns trilinear interpolation of a point (x, y, z) in a unit box
    real(prec), intent(in) :: xyz(3)  ! The point in a unit box [0,1] x [0,1] x [0,1]
    ! V(i, j, k) ... value at vertex i, j, k of the box, where i,j,k = 1,2:
    real(prec), intent(in) :: V(:, :, :)
    real(prec) :: x, y, z
    x = xyz(1)
    y = xyz(2)
    z = xyz(3)
    Vxyz = &
        V(1, 1, 1) * (1-x)*(1-y)*(1-z) + &
        V(2, 1, 1) *   x  *(1-y)*(1-z) + &
        V(1, 2, 1) * (1-x)*  y  *(1-z) + &
        V(1, 1, 2) * (1-x)*(1-y)*  z   + &
        V(1, 2, 2) * (1-x)*  y  *  z   + &
        V(2, 1, 2) *   x  *(1-y)*  z   + &
        V(2, 2, 1) *   x  *  y  *(1-z) + &
        V(2, 2, 2) *   x  *  y  *  z
end function
! *******************************************************


end module searchnew


